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ANALYSIS OF HARMONIC SPLINE GRAVITY MODELS FOR VENUS AND MARS 


INTRODUCTION 

Gravity anomalies are important for exploring the interior structure and 
mass distribution within a planet, and the mechanical properties of its outer 
crust in supporting observed topography. The patterns of gravity anomalies 
were initially analyzed by using sets of spherical harmonic coefficients 
calculated from orbit perturbation data. Examples for Venus include Ananda et 
al. (1980); Bills (1982); Williams et al. (1983) and Mottinger et al. (1983). 
These solutions include information on the long wavelength portion of the 
gravity spectrum only. 

More detailed information on planetary gravity fields is contained in the 
Doppler residual data from radio tracking of spacecraft. The acceleration 
(gravity) anomaly estimates determined from Doppler residuals are for 
components of the gravity field directed along the spacecraft-earth line, 
providing Line-of-sight (LOS) acceleration estimates. A method for estimating 
radial gravity anomalies from LOS acceleration data, their interpolation, and 
use of iteration for improved radial estimates was presented by Bowin (1983). 

A preliminary radial gravity anomaly map of Venus at spacecraft altitude was 
presented, and mass anomalies associated with topographic features were shown 
to be the most significant contributors to the gravity field of Venus. That 
study was supported by NASA Grant NAG 2-91. 

In search of an appropriate means to invert the first radial gravity 
estimate to geoid anomalies for the proposed iteration method (Bowin, 1983), 
we encountered the technique of harmonic splines which allows direct 
estimation of the complete gravity field (geoid, gravity, and gravity 
gradients) everywhere over the planet's surface from the LOS vector data. A 
paper by Bowin, Abers, smd Shure (1985) presented the harmonic spline results 
for Venus as a series of maps at spacecraft and constant altitudes. Global 
(except for polar regions) and local relations of gravity to topography were 
described. That study was supported by NASA Contract NASW-3770. 

HISTORY OF EFFORT 

At the end of the contract period for NASA Contract NASW-3770 we realized 
that corrections were possible for the errors in the zero reference level in 
the LOS values produced by the least squares filtering procedure. We 
submitted (September, 1984) a proposal to NASA to use the harmonic spline 
gravity models as input to the JPL ORBSIM program, to compute the expected 
distorted LOS vectors that would result if the harmonic spline model were the 
actual field. The ORBSIM program allows computed values along an orbit path 
to be integrated along that orbit path to produce distortion of the model's 
gravity field like that distorting the LOS vector values. This procedure is 
described by 
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Phillips et al. (1978) for previous analyses using point and disk mass 
models. The difference between the simulated LOS vectors and the observed LOS 
vectors would then be used to calculate modified LOS vectors. A harmonic 
spline model would then be computed again from the modified LOS vectors. The 
process is then repeated: use ORBSIM to compute distorted LOS vectors from the 
new harmonic spline and compare with the original LOS vectors. Using the 
difference, modify the observed LOS vectors and recompute another harmonic 
spline model from the modified LOS values. The original LOS vectors, of 
course, are retained unblemished for comparison with ORBSIM results of each 
new model. This would be iterated until satisfactory convergence is achieved; 
the resulting harmonic spline model should then represent a close 
approximation to the actual anomalous gravity field of the planet, and its 
distorted product should match the LOS vectors. A variant of this procedure, 
which presumably may permit a faster convergence, is to compute harmonic 
spline solutions for the difference vectors (between simulated and observed 
LOS vectors), and to use those models to estimate the distortion. New 
calculated distorted LOS vectors would then be computed from the combination 
of the distortion harmonic spline model and the original harmonic spline 
model. This process would then be iterated to satisfactory convergence. 

That proposal, however, was declined. The proposal was resubmitted for 
the March 1985 Panel Review. Following discussions in January and February 
1985 with the Discipline Scientist and a member of the former Panel it was 
implied that the prior panel's main criticism was that the ORBSIM procedure 
might not converge because of errors in dynamical parameters used in the 
original processing. I contended, however, that these errors would apparently 
mostly affect the long wavelength values. We discussed this possibility, and 
it seemed to me that, if this proves to be a real problem, spherical harmonic 
model solutions could help damp these errors and their effects. The panel 
member agreed that the panel's "axe" had swung too far in this case. He was 
also surprised to learn that this new methodology was only conceived about two 
weeks before submitting the proposal to meet the deadline. Therefore he could 
then understand why demonstration test cases were not presented in the 
original proposal, the lack of which apparently helped bring the panel to its 
position. Also, at that time, our NASA planetary funds were exhausted: 
payment of the page charges for publication of our paper (Bowin et al., 1985) 
on the gravity field of Venus published in the J6R Proceedings of the 
Fifteenth Lunar and Planetary Science Conference had to be covered by another 
unrelated grant. Because a three-year proposal without proof of the method 
might again be hard for a panel to concur with, I indicated my willingness to 
accept an initial amount that would allow me to test and demonstrate the 
method before NASA need commit to a three-year project. 

On 1 July 1985, I received notification that within 60 days a negotiator 
would negotiate a contract at the $37,696 level for the period May 15, 1985, 
through November 30, 1985. These funds were awarded in late July 1985, and 
the investigation began. Appendix A reproduces an abstract prepared in 
November 1985 which further discusses the 'zero reference level' error and 
contains two illustrations to help the reader better understand the problem. 

A status report of our efforts to October 21, 1985 is best summarized in 
the progress report reproduced on the following three pages. It outlines our 
reasons for obtaining a new minicomputer system to meet the computational 
needs of this and other investigations. 
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PROGRESS REPORT 

NASA Research Grant No. NAGW-792 Entitled 
ANALYSIS OF HARMONIC SPLINE GRAVITY MODELS FOR VENUS AND MARS 
Carl Bowin, Principal Investigator 


An initial amount of $37,969 was awarded in late July 1985 to aid a demon- 
stration test of our proposed new methodology utilizing harmonic splines for 
determining the true gravity field from Line-Of -Sight (LOS) acceleration data 
from planetary spacecraft missions. As is well known, the LOS data incorporate 
errors in the zero reference level that appear to be inherent in the processing 
procedure used to obtain the LOS vectors. Our proposed method offers a solution 
to this problem. 

Considerable progress has been made since the start of this grant, primarily 
in the conversion of programs and requisite data libraries from the Institu- 
tion's VAX 11/780 to a Ridge computer system in our own laboratory. We 
obtained and installed the Ridge to overcome the computational dilemma we had 
experienced 

for the previous two years (i.e., the Institution's VAX charges were rising at 
about 10% per year, the funding for our budgets had remained nearly constant, 
but our computing needs had risen exponentially). The high costs of the VAX 
were particularly traumatic one month when I learned that it cost $2,000 just 
to be logged onto the VAX 8 hours a day, 5 days a week, for one month. Because 
of these problems, I had been exploring the new generation of 32-bit super 
microcomputers and minicomputers. From that investigation we found that our 
needs, the modest equipment budget I had been able to accumulate from industry 
sources over the prior two years, and the design specifications for the Ridge 
computer were very compatible. Thus, in January 1985, we acquired a Ridge 32C 
personal mainframe computer with Reduced Instruction Set Computer (RISC) archi- 
tecture. This computer provides us with a powerful computing resource as re- 
quired for the proposed investigation. 

Our conversion from VAX/VMS to the UNIX operating system, and from VAX-11 
extended Fortran to ANSI Standard Fortran 77 has been a learning exercise but 
tractable. It has taken time and a willingness to work around problems, such 
as a lack of byte variables in standard Fortran 77. And then there are the 
rewards, such as finding that our optimum estimation runs, which took about 6 
seconds per grid point on the VAX 11/780, now take about 4.0 seconds or less on 
the Ridge 32C, and of course, the delight in now being able to undertake 50 to 
100 hour 'experimental' runs at no additional cost. However, we have needed to 
convert several hundred programs and subroutines which is not only time consum- 
ing, but salary consuming. Many of these routines are for general utilities 
such as graph, profile, and map plotting functions, and database manipulation 
routines needed by this investigation. Thus, this progress report is also a 
request for supplemental funding to aid the completion of our demonstration 
test of the new harmonic spline methodology. The VAX version of ORBSIM and 
LUNRES PROGRAMS provided by R. Wimberly of JPL remain on the Institutiton's VAX 
11/780, however, as described later, ORBSIM and LUNRES are now known not to be 
appropriate. 
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The Ridge computer system came without output display devices, so we have 
had to add them in order to have analysis capabilities. One particularly useful 
display device is a 36" wide digital plotter which we are leasing. This plot- 
ting capability is important for testing and evaluating the results of the pro- 
posed harmonic spline modeling methodology for the correction of LOS vector 
data. 

We presently have created on the Ridge system a library of Pioneer Venus 
LOS vectors now organized by revolution number (rev) and by time within each 
rev. This library was made by resorting the LOS library organized by latitude 
and longitude cells used by our prior harmonic spline analysis (Bowin, Abers, 
and Shure, Gravity Field of Venus at Constant Altitude and Comparison with 
Earth, JGR, vol. 90, Supplement, p. C757-C770, 1985). This procedure took 
advantage of the editing already accomplished for that study. We have also 
created a library of basis points for the intended rev-by-rev processing. The 
basis points will be the nodes used in the harmonic spline solutions for each 
rev. Another file identifies the sets of revs to be used for the harmonic 
spline processing. We have started the conversion of harmonic spline coding 
from the VAX-11 extended Fortran to ANSI Standard Fortran 77. The VAX version 
utilized keyed-access data storage, which is not available in standard Fortran. 

We also have examined the functionality of the programs ORBSIM and LUNRES 
and find that neither is appropriate as they now stand for our proposed method- 
ology. Basically, ORBSIM calculates LOS accelerations along given trajectories 
over an input configuration of point or disk masses for comparison with the 
observed LOS accelerations. LUNRES calculates Doppler residuals from the raw 
observed Doppler signal, and then LOS accelerations by diferentiating a smoothed 
cubic spline fit to the Doppler residuals along individual revs. LUNRES is not 
prepared to reproccess a product of LUNRES, ie. LOS accelerations, nor even 
their integral as Doppler residuals. Thus, we find we must now prepare two 
programs. One program will calculate Doppler residuals by integrating a Rev 
time series of LOS accelerations, for the line-of-sight vector, calculated from 
harmonic spline coefficients. The other program will fit a smoothed cubic 
spline to those Doppler residuals, and calculate LOS accelerations from them as 
originally done by LUNRES. We plan to adopt the appropriate subroutines and 
functions from LUNRES for the second program. We also find that the Pioneer 
Venus LOS data we received from S. Ritke of JPL in 1981 has only time tags to 
minute resolution. Hence, time in seconds must be added to our PVO data in 
order to accomplish the integration to Doppler residuals. For the demonstra- 
tion test we expect to hand edit them into the data for the few test revs over 
Beta Regio, and these values will be read from the microfilm that accompanied 
the LOS data tape received from JPL. In a recent telcon, Bill Sjogren of JPL 
has offered to provide a new magnetic tape with LOS time tags to second resolu- 
tion to avoid transcription errors and the enormous effort such hand editing 
would take. This will surely be warranted for further processing following the 
demonstration test success, but to speed our testing we may edit in labels for 
time tags to seconds by referencing the microfilm information. 

The harmonic spline model objectively provides the smoothest model that 
matches input vector LOS data. From these model coefficients we can calculate 
vector data at any selected location in space. It serves a similar role to the 
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point and disk masses used in ORBSIM analyses, but is not dependent upon a 
person for selection of the location and magnitude of point and disk masses, a 
very significant advantage. Our plan therefore is to use harmonic spline model 
coefficients as a memory. In the first iteration, that memory matches the LOS 
accelerations along the orbit trajectory, but upon a lower trajectory (ie., 
downward continuation, errors in the LOS acceleration values become strikingly 
evident. For example see Figure 6 in our paper cited above. Thus the calcu- 
lated LOS accelerations along such a lower and constant altitude can be inte- 
grated to Doppler residuals, and those residuals least-squares fit by a smooth 
cubic spline, which when differentiated will yield calculated LOS accelerations 
with a new zero mean refernce level. These new LOS accelerations will be used 
to determine another harmonic spline model which can then be used to calculate 
it's LOS accelerations along the spacecraft orbit. Those accelerations would be 
integrated to Doppler residuals which would again be smooth spline fit and 
differentiated to calculated LOS accelerations. This process can be repeated 
as desired. Alternately, LOS differences can be used in the computation of 
harmonic spline model coefficients and those coefficients added or subtracted 
to the original set to achieve faster convergence. The object being to obtain 
a harmonic spline model whose LOS accelerations along the rev trajectory when 
integrated to Doppler residuals and then cubic spline fit and differentiated 
and zero meaned by least squares match the original LUNRES determined LOS 
accelerations. At that stage the harmonic spline model provides an estimate of 
the true planetary gravity field beneath that rev path to the data resolution. 
Such rev harmonic spline models then can be used to construct grids of gravity 
values at constant altitudes over the planet surface. Such grids are vital for 
admittance studies of planetary structures by gravity/ topography transfer 
function methods, and offer the strongest means by which to examine the 
internal structure of another planet. 


HISTORY OF EFFORT CONTINUED 


The preceding PROGRESS REPORT was also a SUPPLEMENTAL FUNDING REQUEST 
UNDER NASA Grant Award No. NAGW-792 (for $27,275) for one month's funding 
support. The Lunar and Planetary Review Panel met in early November 1985 and 
recommended that no additional funding should be considered until the pilot 
project is completed. My letter to the Discipline Scientist of November 22, 
1985 indicated that I could understand the reasons for such a decision, and I 
could also understand (and feel) the reasons for our needing funding. To have 
adequate computing resources for the Venus and Mars gravity task, and other 
projects, I had to obtain our own high-performance computing resource. The 
conversion of software to it has taken more time that originally thought (a 
common ailment). We also had to develop our own routines, for example for 
matrix inversion, which formerly were provided as commercial IMSL routines on 
the Institution's VAX. IMSL was willing to allow us to port the matrix 
inversion routine to the Ridge computer but at a cost of about $3,000. Our 
efforts were also required to add time in seconds to the LOS vector data for 
the PV data from JPL so that it can be integrated back to Doppler data. Thus, 
that letter requested a no-cost extension to 30 June 1986 for our grant which 
was to have terminated November 30, 1985. 

Following conversion of our PV LOS library to rev-by-rev with time tags, 
and conversion of our harmonic spline programs from the VAX 11/780 to the 
Ridge 32C computer and its modification for the processing of individual and 
adjacent revs, development of a node database for such processing, testing 
began. Unfortunately, we could not get the ensuing matrix to invert without 
error. On January 7, 1986 I requested a no-cost extension to 31 December 1986 
for Research Grant No. NAGW-792, for which approximately $6,000 remained. 

We solved a problem with our matrix inversion routine that improved 
inversion of the data matrices used in our Optimum Estimation program for 
global earth studies. We hopefully thought that change would also benefit the 
harmonic spline processing of PV revs, but it did not. None of our experi- 
ments through 31 December 1986 obtained a successful matrix inversion for a 
single rev supplemented by data for the two adjacent revs. Without this 
important inversion, testing could not proceed, and our funds were exhausted. 

We tried to the best of our ability and resources to prove our new 
suggested methodology for determining the true detailed gravity field of Venus 
and Mars from LOS data. I have additional ideas in mind that may help to 
resolve this dilemma and thereby obtain an inversion for the needed harmonic 
spline model. I will try them if we have the opportunity to self-fund the 
effort in the future. The recent increased difficulty in obtaining federal 
research support, the soft-money status of research at W.H.O.I., and the lack 
of encouragement from NASA, however, makes this prospect dim. 
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Appendix A 


ANALYSIS OF HARMONIC SPLINE GRAVITY MODELS FOR VENUS AND MARS 
Carl Bowin, Woods Hole Oceanographic Institution, Woods Hole, MA 

The use of harmonic spline methodology to resolve planetary gravity 
fields from the randomly oriented residual line-of-sight (LOS) vectors that 
result from processing Doppler tracking data was first demonstrated in our 
work published earlier this year (Bowin et al. , 1985). The LOS vectors 
are components of the planetary gravity field along that direction that 
happened to point toward the earth at the time of measurement. The LOS 
acceleration vectors are computed by differentiating the residual Doppler 
signal that remains after the Doppler contribution from the central body 
attraction, along the computed orbit trajectory is subtracted from the 
observed Doppler data. Least-squares filtering procedures, applied to the 
Doppler data to obtain the LOS results reduce the true absolute amplitudes 
of the acceleration vectors (Gottlieb, 1970), perhaps by as much as 40 
percent. Apparently that reduction of a maximum of a positive anomaly is 
accompanied by the pumping of energy into lows that flank a gravity high. 

I view the reduction of the positive amplitude as most likely, indi- 
cating an uncertain location of the zero reference level. For example, 
imagine a volcano or ridge upon a planet with a very strong crust. In 
this case the topographic load would be completely uncompensated, and the 
free-air gravity anomalies across the structure would be entirely posi- 
tive, as shown by the solid curve, in reference to the dotted zero line in 
Figure 1. An integral of the gravity acceleration curve would yield 
anomalous Doppler frequency variations. Upon a least-squares fitting of a 
smooth curve through such Doppler data, and its subsequent differentia- 
tion, an anomalous gravity curve shown by deviations from the dashed line 
(which represents an inferred erroneous zero reference level) might 
result. The resulting gravity profile might then have a reduced positive 
anomaly with flanking lows. In such a situation, the integral of the 
positive and negative gravity anomalies along the orbit path would tend to 
approach closer to zero than in the true condition. 


+ 



Figure 1. Hypothetical gravity anomaly. See text for 
discussion. 
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The above hypothesis is supported by the occurrence of gravity lows 
flanking positive gravity anomalies associated with topographic highs on 
Venus and Mars. Where the topographic high is a nearly circular feature, 
the lows symetrically straddle the positive along the orbit tracks. This 
is true even though the spacecraft altitude at the two symmetric lows may 
be quite different. See Figure 2 for the example of this relationship in 
the orbit passing above the Crest of Beta Regio on Venus. This observa- 
tion is seen also to be true in maps of LOS accelerations at spacecraft 
altitude for Venus published by Sjogren et al. (1980,1981), Bowin (1983) 
and Bowin et al. (1983) for example. In fact it was these relations that 
led to the above hypothesis. 


Tracking data and spacecraft 
ephermeris have been used to com- 
pute the coefficients for spherical 
harmonic expansions of the poten- 
tial field of planets. For example, 
see Bills (1982) and Mottinger et 
al. (1983) for Venus. These global 
solutions only have low resolution, 
however. Point masses and disk 
models over a grid have been used 
to compute Doppler residuals along 
orbit trajectories by the program 
ORBSIM and results compared with 
the observed LOS residual accelera- 
tions. This process has allowed 
the construction of radial gravity 
anomaly maps from the point mass or 
disk models. Depending on where or 
how the zero reference errors arise 
during processing, they may remain 
in the point and disk mass models 
since they match to the LOS data. 

That the lows flanking the gra- 
vity highs are in large part arti- 
facts from the LOS processing is 
clearly documented by the harmonic 
spline results (Bowin et al. , 

1985). The power of the harmonic 
spline methodology is that it ob- 
jectively solves for the smooth- 
est field that will match the 
observed data, and that the data 
may be of random components of 
the field, and at random alti- 
tudes. Then from the harmonic 
spline coefficients at the basis 
points used for the computation, 
any component of the gravity 




S DISTANCE (km) N | 

Fig. 6. Comparison of profiles of ra- 
dial gravity anomaly at spacecraf t, 500 
km,& 300km altitudes across Beta Regio 
in a direction parallel to the orbit 
paths. The line types in the height 
plot (A) provide identifications for 
the lines describing the gravity pro- 
files (B) (Bowin et al, 1985, Fig. 6). 
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field at any location and altitude may be determined. In addition to the 
gravity profile at spacecraft altitude in Fig. 2, we show the curves at 
constant altitudes of 300 and 500 km. Note how the low north of Beta 
Regio becomes greatly exaggerated upon increasing downward continuation. 

In the 300 km map (Bowin et al. , 1985, Fig. 5) that low is seen to greatly 
exceed the absolute magnitude of all other anomalies on Venus, and thus it 
is obviously in error. Because downward continuation is a rigorous proce- 
dure, that error signifies that its value at spacecraft altitude is too 
negative. Note, also, that the low south of Beta Regio is upward continued 
to 300 and 500 km and thus is attenuated at those constant altitudes. 
Through the use of harmonic spline models, in an iterative manner, we are 
now attempting a procedure by which a harmonic spline model can be obtained 
that when appropriately filtered will match the observed LOS data along 
the spacecraft trajectories, and when unfiltered j" will represent the true 
anomalous gravity field of the planet. At that stage, we can construct 
grids of gravity values at constant altitudes over the planet surface. 

Such grids are vital for admittance studies of planetary structures by 
gravity/topography transfer function methods, and offer the strongest 
means by which to examine the internal structure of another planet. 

For this computer intensive study, as well as for other computing 
needs, we have acquired a Ridge personal mainframe computer with RISC arch- 
itecture. This system provides us with an essentially unlimited computer 
resource. The 'noise' in the harmonic spline results of Fig. 2 along the 
spacecraft trajectory result from differences in solutions at boundaries of 
latitude and longitude cells previously used to divide the large LOS data 
set into manageable sizes. For the present study we have created on the 
Ridge system a library of Pioneer Venus LOS vectors organized by revolu- 
tion number (rev) and by time within each rev. Harmonic spline models will 
be computed for entire revs, thus eliminating the noise described above. 
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